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Abstract 

The quasi-periodic doubling cascade is shown to occur in the transition from regular 
to weakly turbulent behaviour in simulations of incompressible Navier-Stokes flow 
on a three-periodic domain. Special symmetries are imposed on the flow field in 
order to reduce the computational effort. Thus we can apply tools from dynamical 
systems theory such as continuation of periodic orbits and computation of Lyapunov 
exponents. We propose a model ODE for the quasi-period doubling cascade which, 
in a limit of a perturbation parameter to zero, avoids resonance related problems. 
The cascade we observe in the simulations is then compared to the perturbed case, 
in which resonances complicate the bifurcation scenario. In particular, we compare 
the frequency spectrum and the Lyapunov exponents. The perturbed model ODE 
is shown to be in good agreement with the simulations of weak turbulence. The 
scaling of the observed cascade is shown to resemble the unperturbed case, which 
is directly related to the well known doubling cascade of periodic orbits. 
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1 Introduction 

In the absence of boundaries the incompressible Navier-Stokes equations are 
symmetric under rotations, translations and reflections. If we impose a sub- 
group of symmetries on the solutions we reduce the number of degrees of 
freedom in simulations of the flow. In the early days of numerical simulation 
such reduction by symmetry was used to probe into turbulent flow, revealing 
structures such as the Taylor-Green vortex pQ. Kida [2] put forward what is 
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probably the maximal reduction still allowing for turbulent flow. The corre- 
sponding flow is called high symmetric and was used to study flow statistics 
at moderate to high Reynolds number [3"Pfo] . 

As computers have grown considerably since then, those results may be repro- 
duced today in simulations without any special symmetry. However, reduction 
by symmetry can still be a useful approach when going to ever higher Reynolds 
number, requiring higher resolution of the models, and when applying tools 
of dynamical systems theory to turbulence. 

Recently considerable effort has been made to apply continuation and bifur- 
cation analysis to fluid dynamical problems, see e.g. [6]Tf] . If this could be 
done at realistic resolution the results could prove a valuable complement to 
statistical analysis of direct numerical simulation. At the moment the maxi- 
mum number of degrees of freedom tackled successfully is of order 10 4 . Given 
that the simulation of turbulent flow at moderate Reynolds number requires 
a number of degrees of freedom in the order of 10 6 , reduction by symmetry 
can be successfully applied. 

Here, we examine the transition from regular to weakly turbulent motion in 
high symmetric flow. Exploiting the divergence free condition in addition to 
the symmetry we gain a factor of about 300 with respect to general, non sym- 
metric flow and a factor of 3/2 with respect to earlier work on high symmetric 
flow in terms of the number of degrees of freedom. Thus we can analyse the 
transition using continuation of periodic orbits and computation of Lyapunov 
exponents, using the direct method for integration rather then the pseudo- 
spectral method. 

In addition to the Ruelle-Takens scenario, reported on previously [8J, we find 
that a cascade of quasi-periodic doubling bifurcations, otherwise knows as 
torus doubling bifurcations, occurs in the transition to weak turbulence. To 
the author's best knowledge this is the first time this bifurcation sequence 
is found in Navier-Stokes flow at realistic truncation level. In order to relate 
the numerical results to theoretical predictions we propose a model Ordinary 
Differential Equation (ODE) for this cascade. In contrast to model equations 
studied previously PTTU] these equations, based on a fixed time smooth sus- 
pension of the Henon map, display a complete cascade with a fixed frequency 
ratio of the bifurcating tori in the limit of a perturbation parameter to zero. 
If we introduce a slight perturbation in the model ODE the bifurcation cas- 
cade is interrupted, as it is in the simulations of weak turbulence. Also, the 
frequency spectrum at the onset of chaos and the dependence of the leading 
Lyapunov exponents on the bifurcation parameter in the simulations is shown 
to be very similar those of the perturbed model ODE. 

To within the accuracy of the numerics, the scaling of the cascade is shown 
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to agree with that of the unperturbed model ODE, which is known theoreti- 
cally. When looking at the leading Lyapunov exponents a difference becomes 
apparent: in the unperturbed model ODE the largest two nonzero exponents 
become equal on an open interval in parameter space in between successive 
doubling bifurcation points, whereas in the perturbed model ODE, and also 
in the flow simulations, they remain separate. This difference might by caused 
by the crossing of resonance tongues or other codimension one phenomena in 
which the stable torus loses normal hyperbolicity. A thorough study of the 
model ODE with finite perturbation might yield an explanation. 



2 The vorticity equation for high symmetric flow 



Consider an incompressible fluid in a periodic box < x\,X2,xs < 2ti. In 
terms of the Fourier representation of velocity and vorticity, 
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we have 



^a)j(k) = etjkkjki v k v t - zyfc 2 d>i(k) (2) 

hui = (3) 
£>i(k) = -ei jk kjV k (k) (4) 

where v is the kinematic viscosity, e^j^ is the permutation symbol and the 
tilde denotes the Fourier transform. In terms of the standard norm energy 
and enstrophy are given by 



e = -\\A 2 Q = ^\\ 2 (5) 

respectively. Now consider the following discrete symmetry operations: Si, 
reflections in the planes Vi given by Xi = ir and R4, rotations over ir/2 radiants 
about the lines U : Xj = tt/2 for j 7^ i. If the flow is invariant under these 
operations only one out of three components of vorticity in a volume fraction 
1/4 3 needs to be computed to determine the flow on the periodic domain. 
Figure [T] gives an impression of the symmetries. We call a flow invariant under 
Si and B4, first described by Kida [2], high symmetric. 

Symmetry operations Si and B4 introduce linear relations between the Fourier 
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Fig. 1. The structure of high symmetry. F represents the state on a face of a box 
with edges of length tt/2. Reproduced from [2]. 

components of vorticity. First of all we have 

fc 2, h) = U2(h, h, h) = h, h) (6) 



so that we may consider only one component. This scalar function is even or 
odd in its arguments: 



wi(fci, k 2 , k 3 ) =u 1 (-k 1 , k 2 , h) = 

-u)i(h, -k 2 , h) = -&i(ki, k 2 , -k 3 ) 

and finally we have 



(7) 



&i{ki, k 2 , k 3 ) = < 



— Oi(ki, k 3 , k 2 ) for ki and k 2 and k% even, 
&i(ki, k 3 , k 2 ) for k\ and k 2 and k 3 odd, 
otherwise. 



We consider a cubic truncation, i.e. I&i^l < N. Relations ([ME]) reduce the 
number of independent Fourier modes of vorticity by a factor of 192 in the 
leading order, that is N 3 . This reduction was exploited by Kida et al. [3|4] to 
investigate scaling laws at moderate to high Reynolds number. In the present 
work we exploit the divergence free condition for vorticity to further reduce 
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the number of modes. With the aid of equation (JS)) it reads 
hcj^h, k 2 , k 3 ) + k 2 &i(k 2 , k 3 , h) + faoj^fa, k h k 2 ) = 



(9) 



Taking maximal advantage of relations ([Ml]) we consider only Fourier compo- 
nents of U\ in the fundamental domain {k G Z 3 |A; 3 > k 2 , k 3 > k\, ki > 0, k 2 > 
0, k 3 < N}. The number of independent modes is reduced by a factor of 288 
in the leading order. 

These components satisfy the following equation 

^wi(/ci, k 2 , fa) = k 2 k 3 (S(k 3 , fci, k 2 ) - S(k 2 , k 3 , k{)) 

+k 1 k 2 f(k 2 , k 3 , fci) - kzkxffa, fei, k 2 ) 

+ (k 2 2 - kDfih, k 2 , k 3 ) - vk 2 u x {k Xl k 2 , k 3 ) (10) 

where S and T are the Fourier transforms of 



S(x 1 ,x 2 ,x 3 ) =v 1 (x 1 ,x 2 ,x 3 ) 2 

T(x 1 ,x 2 ,x 3 )=v 1 (x 2 ,x 3 ,x 1 )v 1 (x 3 ,x 1 ,x 2 ) (11) 

and 

k 2 v x = k 2 tii{k 3 , ki, k 2 ) - k 3 ^i{k 2 , h, k x ) (12) 



Energy in supplied by fixing the low order odd mode o)i(l, 1, 3) = —3/8. Thus 
we obtain a family of dynamical systems with one parameter, u, and a number 
of degrees of freedom given by 

, l(f) 3 + l(f) 2 -If - 1 if AMseven / x 

n(N) = I 3 V 2 > o 2 V 2 > 62 13 
2 (n-i^ _|_ i fiLzl) _ I (m^l) - 1 if AT is odd. 
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3 Numerical considerations 



In performing time integrations we avoid the use of a pseudo-spectral method, 
commonly employed for three dimensional simulations of Navier-Stokes. Up to 
a truncation level of iV = 25 (n = 1365) the direct method is not exceedingly 
slow due to the reduction of the number of degrees of freedom described above 
and yields easy access to the Jacobian for integration of the variational equa- 
tions. Also, the direct integration code can easily be run in parallel, distribut- 
ing the computation of the components of the vector field over any number 
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of processors. For higher truncation levels, at which storage of the nonlinear 
interaction coefficients requires huge memory space, it is mandatory to use a 
pseudo spectral method, as described in [3J, even though it does not exploit 
the divergence free condition . For time integration a seventh to eight order 
Runge-Kutta-Felbergh scheme with adaptive step size is employed. As many 
of the results presented here are based on rather long time integrations it 
important to keep the error tolerance low. We have checked energy conserva- 
tion at zero forcing and viscosity using the high order and the fourth order 
Runge-Kutta schemes. For realistic 0(1) levels of the energy, an integration 
time At = 10 3 and a fixed error tolerance 8E = 10~ 9 the step size required 
for the fourth order scheme is about ten times smaller then the average step 
size using the high order scheme. As the high order method needs thirteen 
evaluations of the vector field at each time step, against four for the low order 
scheme, it is more efficient. 

Below simulations are performed with a viscosity in the range 0.005 < v < 0.01 
and the truncation level fixed to N = 15 (n = 300). We computed the energy 
and the enstrophy, as well as Taylor's micro-scale Reynolds number, R\, and 
Kolmogorov's dissipation length scale, 77, defined by 



The dimensionless number r]N indicates if the resolution is high enough in 
numerical simulations. If rjN ~ 1 the truncation error is considered negligible. 
At the fist transition to chaos, the focus of this paper, we have 7]N w 0.8. The 
band average energy spectrum is shown in figure [2]( left). In order to make sure 
that the results presented below do not depend critically on the truncation 
level we repeated some of the computations, in particular the location of the 
quasi-periodic bifurcation points described in section [5j at iV = 21 (n = 814). 
The qualitative behaviour remains the same as the bifurcation points shift 
slightly to lower viscosity. For this truncation level we have r]N rs 1 at the 
first transition point. The band-averaged energy spectrum is shown in figure 
fright). At the small scales an exponential decay is visible, indicating that 
our numerical results are reliable. 

Periodic orbits are continued in the viscosity as fixed points of a Poincare 
map, using the arclength method as described in [IT] . This method is time 
consuming because of the integration of the linearised equations. However, 
running the code in parallel this is no obstacle. 




(14) 
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Fig. 2. Band-averaged energy spectra near the first transitions to chaos, described 
in section Obtained from an integration of 300 time units, log-linear scale in 
normalised units. Left: truncation level N = 15, right: N = 21. 



4 Transitions to chaos 



For large viscosity an equilibrium state is the global attractor in our simu- 
lations of high symmetric flow. At v ps 0.0113 a Hopf bifurcation occurs in 
which a stable periodic orbit is created. To get an impression of the transi- 
tions from periodic to chaotic motion we computed a limit point diagram in 
the range 0.01 > v > 0.005 of the Poincare map on the coordinate plane given 
by u(0, 2,4) = —0.05. In this range the time average micro-scale Reynolds 
number varies from R\ p=i 50 to R\ ~ 27, indicating that the aperiodic be- 
haviour can be classified as weak turbulence, fully developed turbulence sets 
in around v ps 0.004 and R\ ps 60 in high symmetric flow. The limit point 
diagram is shown in figure [H Two parameter ranges with periodic behaviour 
can be seen: one around v = 0.01 and one around v = 0.0068. A continuation 
of these periodic orbits in parameter v is shown in figure HI The branch which 
is stable around v = 0.0068 does not bifurcate from an equilibrium at high 
viscosity. Both branches become unstable in a Neimark-S acker bifurcation. 
Directly beyond these bifurcation points we expect the behaviour to be quasi 
periodic, and indeed invariant circles appear in the Poincare section, figure [3j 
The breakdown of these invariant tori gives rise to chaos. The torus created 
at NS2 displays a quasi-periodic Hopf bifurcation and thus the Ruelle-Takens 
scenario to chaos is followed, as reported in [8]. The torus created at NSi dis- 
plays a quasi-periodic doubling bifurcation which turns out to be the first of 
a cascade, discussed in detail in section |5j 

The fundamental frequencies of the tori created at NSi,2 have a physical inter- 
pretation. The higher frequency, U\ ~ 2.6 is set by the large eddy overturning 
time, estimated as T = l/U ~ 2.3 where U ~ 0.9 is the root mean square 
velocity and I = <\/it 3 /3 is the characteristic length scale given that by sym- 
metries Si, the planes £1,2,3 — tt are impermeable and by symmetries Ri the 
velocity field has a three fold rotational symmetry around the main diagonal 
X\ = x 2 = £3. The lower frequency, w 2 ~ 0.26, corresponds to a modulation 
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of the amplitude of the energy. It is also present in fully developed turbulence 
with high symmetry and can be though of as a retaining time scale of anoma- 
lously high energy. The dynamics on these time scales is illustrated by the 
isosurface of enstrophy shown in figure On the short time scale, T, patches 
of high vorticity are generated on the main diagonal and off the diagonal 
in triples related by the three fold rotational symmetry. These patches con- 
verge to the origin and to the center of the 27r-periodic box, points of special 
symmetry. On the longer time scale the overall amplitude of this process is 
modulated. Note, that the patches of high vorticity do not have the elongated, 
tubular shape typical of developed turbulence. 



5 The quasi-periodic doubling cascade 

In figure [6] six Poincare plots are shown beyond bifurcation point NSi. The 
torus "doubles" at least four times (a-e) before a chaotic attractor appears 
with a shape very similar to the doubled torus. 

Such sequences of torus doubling bifurcations have been observed both nu- 
merically, e.g. in severe truncations of the Navier-Stokes equations [T7] and 
in a periodically driven low-order atmosphere model [T5], and experimentally, 
e.g. in electronic circuits [T5][T6] and have also been studied in connection to 
strange non chaotic attractors [19] . In these examples chaotic behaviour is 
observed after two or three doublings. In our system we find four doublings, 
after which quasi-periodic motion is very hard to distinguish from aperiodic 
motion numerically. 

The observation of incomplete quasi-periodic doubling cascades led to the in- 
troduction of simple models which highlight the similarity to the well known 
periodic doubling cascade [9"][T0] . In these models an infinite cascade can only 
be observed in a limit of the parameters were the vector field can be decom- 
posed into a three dimensional vector field which exhibits a period doubling 
cascade and a constant rotation. Away from this limit the stable torus is more 
fragile after each bifurcation and the cascade is interrupted. 

In order to see why the cascade is interrupted, we need to investigate the 
quasi-periodic doubling bifurcation in detail. Rigorous bifurcation theorems for 
bifurcating invariant tori were formulated by Broer et al [T2iri3] . They include 
the bifurcation shown here, in which a torus loses stability and a new stable 
torus is created with one fundamental frequency half that of the original torus. 
However, the bifurcation theorem is given for a one parameter family of tori 
that satisfy a non resonance condition on the fundamental frequencies. If the 
fundamental frequencies of the torus become resonant normal hyperbolicity is 
lost and persistence under parameter variations is not guaranteed. We have 
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only one parameter in our system and the fundamental frequencies depend 
on it continuously. As we vary the parameter we cut through a Cantor set 
of values at which the torus is normally hyperbolic. Thus, the loss of normal 
hyperbolicity through resonance is a possible mechanism for the interruption 
of the cascade. 

However, the resonances we pass through, i.e. the Arnold' tongues we cross, 
are of fairly high order as we find for the fundamental frequencies at NSi that 
oj i ~ 10u2- Recently an algorithm for the continuation of invariant tori was 
developed that "steps over" such high order resonances without a problem, 
and in fact this algorithm was tested on a system that exhibits a cascade of 
quasi-periodic doublings [H]. In this paper it is pointed out that in the absence 
of a non resonance condition the notion of a quasi-periodic bifurcation point 
itself is unclear, as no normally hyperbolic torus exists in a whole interval in 
parameter space. On one side of this interval we observe the stable single torus 
and on the other side the stable doubled torus. This interval turns out to be 
small compared to the numerical resolution of our simulations, a step size in 
the viscosity of about 5u fh 1CT 7 . 

Apart form resonances there are other mechanism that can interrupt the cas- 
cade. Stagliano et al. [19], for instance, concluded that the cascade was inter- 
rupted when the stable torus loses smoothness and collides with an unstable 
parent torus. In reference [T5] a strikingly complicated, yet incomplete picture 
is sketched of the mechanisms that can interrupt the cascade. 

In order to interpret our numerical results, let us formulate a new model ODE 
for the quasi-periodic doubling cascade. Models studied before [9"][T9"] were 
based on the logistic map and cannot be regarded as the Poincare map of a 
smooth flow. Our starting point is the Henon map. 

Let H s (x, y) be a family of Henon maps on R 2 which exhibits a period doubling 
cascade. Let Sj be the period doubling bifurcation point at which a stable T- 
periodic point is created and let lim^oo gj = s c . Let <f>t(x, y, it; s) be a smooth 
fixed time suspension of l-i s on I 2 x S so that <fii(x, y, 0; s) = (H s (x,y),0). 
The flow (fi t is generated by an autonomous vector field, denoted here by 
fi9 x + hd y + d u . In appendix IA1 explicit expressions are given as formulated 
in [20] . Finally, we add a second periodic variable and a coupling proportional 
to e: 

x = fi(x, y, u; s) + ehi(v) u=l 



with h : S ->■ R 2 , g : R 2 ->■ R and u irrational. For e = the flow on R 2 x»SxS 
is given by ijjt{x, y, u, v; s) = (4>t( x , V, u \ s \ + ut) and we know that 





(15) 
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(1) At s = Si a stable torus with fundamental frequencies 1/2* and u> is 
created. 

(2) The quasi-periodic doubling bifurcations accumulate on zero with the 
same scaling as for the period doubling cascade, thus for the bifurcation 
points we find that Hmj_ ) . QO (sj — Sj_i)/(sj + i — S{) = 6, where 5 = 4.669 . . ., 
the Feigenbaum constant. 

(3) The Lyapunov exponents of the attracting torus are determined by the 
Floquet multipliers of the corresponding periodic points of H. s . Therefore 
the leading nonzero exponent should be equal to the second exponent on 
a open interval contained in each interval (sj, s i+1 ). 

For finite e the second fundamental frequency will depend on s continuously 
and resonance points will be passed through on the way to the accumulation 
point s = s c . 

In the following we compare the model ODE to the simulations of weak tur- 
bulence. The damping parameter of the Henon map is fixed to j = 1/4. The 
unperturbed model then displays a doubling cascade with s± = 1.1718..., 
S2 = 1.7031 . . ., S3 = 1.8144 . . . and s c = 1.84525 .... In the perturbed system 
we fix the second frequency to the golden mean, i.e. u = (\/5+ l)/2 and 



h\{v) = sin27ru 
h 2 (v)=0 

9(?,y)=y (16) 

The perturbation parameter is fixed to e = 0.001. At this perturbation strength 
we observe four doubling bifurcations before chaos sets in. In the Poincare in- 
tersection plane given by v = 0, projected on the x and y variables, the 
bifurcation sequence looks very similar to the one shown in figure El 

We have computed the frequency spectrum near the onset of chaos and the 
(partial) Lyapunov spectrum both for the weak turbulence and for the model 
ODE. For computation of Lyapunov exponents of the high-symmetric flow we 
used finite differencing rather than integration of the full linearised system. 
Integration was done until the zero exponent associated with the direction 
transversal to the flow and tangent to the torus had converged up to 10~ 4 . 

Turning first to the frequency spectrum we observe that after three doublings 
bifurcations peaks at frequencies U2/2, Co>2/4 and U2/8 show up, as well as 
combination peaks around u\. Figure [7] shows the spectrum around the fun- 
damental frequencies for the simulations of weak turbulence. A rather similar 
figure for the model ODE is omitted here. In the case of a period doubling cas- 
cade the amplitude ratio \i between consecutive peaks u/2 l can be predicted 
on basis of the scaling constants [21] . We attempted to measure this ratio for 
the spectrum shown in figure [7] but only the first three peaks can be distilled 
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from these data. These scale roughly by a factor of fi = 4, not far from the 
theoretical prediction for the period doubling case /i = 4.648 .... The same 
estimate is obtained for the model ODE. In the case of Poincare maps, rather 
then discrete maps, similar deviations from the theory are found. 

In figures [Deleft) and deleft) the leading Lyapunov exponents are shown as 
a function of the bifurcation parameter, i.e. the viscosity v for the high- 
symmetric flow and the parameter s of the Henon map for the model ODE. 
From these data we can estimate the locus of the bifurcation points, stressing 
once more that this notion is ill-defined on a scale smaller then our resolu- 
tion in parameter space. Thus we find consecutive estimates for the scaling 
constant 5 as 5 { = {pi - zVi)/Oi+i - Vi) or 5 { = (s s - Sj_i)/(s i+1 - Sj): 

5i = 4.64 ± 0.2 $2 = 4.88 ± 0.5 for high-symmetric flow and 
5 1 = 4.75 ± 0.05 5 2 = 4.63 ± 0.2 for the model ODE. 

Again, this is compatible with the theory for period doubling cascades, albeit 
a rough estimate. 

In figure fright) enlargements around the points D x and D 2 are shown for the 
high-symmetric flow simulations. After the first doubling the leading non zero 
exponents are equal on small interval of order 8v = 10~ 6 . After the second 
doubling the leading exponents neither cross nor become equal. This result was 
confirmed by integration of the full linearised equations at several points in the 
parameter range where the leading exponents attain their minimal distance. 
The parameter range between the third and fourth doublings is too small to 
see whether the exponents meet. In this range we record the passage through 
a high order resonance which does not destroy the torus. 

In figure fright) an enlargement is shown of the Lyapunov exponents of the 
model ODE around the doubling points D3 and D4. For comparison, the ex- 
ponents of the unperturbed model, computed from the Henon map itself, 
have been drawn with thin lines. Again the Lyapunov exponents do not be- 
come equal in between doubling points. The step size has been taken as small 
as 5s = 1 • 10" 5 . Apparently, the slightly perturbed model ODE yields the 
same behaviour of the leading Lyapunov exponents as the simulations of high- 
symmetric flow. 

For now we have no explanation for the absence of crossing or equal leading 
Lyapunov exponents in between doubling points. As noted before, the invari- 
ant tori are structurally stable only on a fractal domain of the parameter and 
bifurcations might go unnoticed when monitoring Lyapunov exponents with 
finite step size for the parameter. A more thorough understanding could be 
obtained by a continuation of the invariant torus itself and by further study 
of the model ODE (fT5l) . which is the topic of future research. 
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6 Conclusion 



We have investigated the transition from stationary to disordered behaviour 
in simulations of high symmetric flow. In combination with the symmetry 
the divergence free condition was exploited to reduce the number of degrees 
of freedom by a factor of about 300 with respect to simulations of general 
periodic Navier-Stokes flow, a gain of 30% compared to earlier work [2]- [5], [8]. 
This allowed us to investigate the bifurcation scenario in detail and at realistic 
truncation level. 

Along with the Ruelle-Takens scenario, the quasi-periodic doubling cascade 
occurs as a route to weak turbulence. By means of Poincare sections, power 
spectra and Lyapunov exponents we have shown that this cascade bears close 
resemblance to the well known doubling cascade for periodic orbits. We observe 
four doublings after which quasi-periodic behaviour is very hard to distinguish 
from chaotic behaviour numerically. In previous work on the quasi-periodic 
doubling cascade, both numerical and experimental, chaos is reported to set 
in after two or three doublings jTo]ll7|16] . The fact that we see more doublings 
might be due to the widely separate fundamental frequencies of the bifurcating 
tori in our system, as is the case in recent work by Schilder et al[14j. 

In order to make the notion of a quasi-periodic doubling cascade precise we 
have proposed a model ODE, equations (ITo"]) . such that in the limit of a per- 
turbation parameter e j. the correspondence to the periodic doubling cascade 
is exact. The scaling of the first steps of the observed quasi-periodic doubling 
cascade agrees with the prediction of the unperturbed equations f fT5|) . Looking 
at the leading nonzero Lyapunov exponents, however, a difference between the 
perturbed model ODE and the high-symmetric flow simulations on one side, 
and the unperturbed model ODE on the other side, becomes apparent. In the 
perturbed model ODE and in the simulations the two largest nonzero expo- 
nents remain separate in between doublings, or at least no crossing is found 
with a small but finite step size in the parameter. An explanation might fol- 
low from a detailed study of the model ODE (|T5|) with nonzero perturbation. 
Other interesting questions would be if the scaling laws of the doubling cas- 
cade are influenced by the perturbation and how many steps of the cascade 
can be observed given the ratio of fundamental frequencies and the strength 
of the perturbation. 

The solutions presented here, in particular the bifurcating tori, are solutions of 
the truncated three dimensional Navier-Stokes equations on a periodic domain. 
They might, however, prove to be unstable to asymmetric perturbations. It 
is therefore uncertain whether the quasi-periodic doubling cascade can be 
observed in periodic flow without any symmetries imposed. It does, however, 
occur at low Reynolds number, at which asymmetric perturbations might be 
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damped. Even at such low Reynolds number, the simulation of general three 
dimensional flow, and the bifurcation analysis as presented here, remains a 
formidable task. 
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A An explicit suspension of the Henon map 

Explicit expressions for a vector field that yields a fixed time, smooth suspen- 
sion of the Henon map were given by Mayer-Kress and Haken [20]. Here we 
repeat the essential definitions in a the notation of the model ODE [15j 

The Henon map is given by 
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In order to define the function f(x, y,u; s) we need to introduce the following 
auxiliary functions: 



f ( u ) = (3 - 2u)u 2 

rj( u ) = l+u\nj - (21nj + 3)w 2 + (2 + \nj)u 3 



i 

i 2 + r] 2 



e+v 2 , 



X(x,y)=j u (r)x-£y) 



(A.2) 



where the primes denote differentiation with respect to -u. Suppressing the 
explicit dependence on u in the right hand side we have 



f(x,y,u;s) 




F(X(x,y);s) 



(A.3) 
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, VI chaotic 

, ' V two or three periodic 
/ , — 
' / IV periodic 




Fig. 3. Limit point diagram of high symmetric flow simulations. Iterates of the 
Poincare map, obtained from an integration of length At = 1000, are drawn for each 
parameter value. Transient motion is discarded and the step size of the parameter 
is set to Su = 1 ■ 10~ 4 . The last point of each integration is the initial point at the 
next parameter value, starting at the stable periodic orbit at v = 0.0105. Visible is 
the transition I— >ll from periodic to quasi periodic, II— to chaotic, III— TV back 
to periodic, IV— >Y to quasi periodic and V— >VI to chaotic/turbulent. In region III 
the behaviour alternates between chaotic and two or three periodic while the spatial 
structure of the flow remains simple. 
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Fig. 4. Branches of periodic orbits of high symmetric flow, truncated at N = 15. 
Solid lines denote stable branches and dashed lines unstable branches. Bifurcation 
points are marked SN for saddle-node, NS for Neimark-S acker and H for Hopf. 
The roman number refer to the stable periodic behaviour as seen in the limit point 
diagram, figure [3l 
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Fig. 5. Isosurface of enstrophy during quasi-periodic behaviour corresponding to 
figure [6^,. Edges of length ir, taken from the 2-7r-periodic domain and isosurface level 
set to eight times the spatio-temporal average. The main diagonal is an axis of 
three-fold rotational symmetry. The time step between a and b (c and d) is half the 
short period T\ = 2tt/ 'oj\ and from a to c (b to d) half the long period Ti = 2tt/uj2- 
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Fig. 6. Poincare sections in the plane w(0, 2,4) = —0.05, projected onto the modes 
w(2,2, 4) and cj(3, 1,3). Four quasi-periodic doublings are shown in (a-e) and the 
resulting chaotic attractor in (f). 
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Fig. 7. Frequency spectrum of the energy just beyond bifurcation point D3 around 
the widely separated frequencies LO± t 2- Note the peaks due to the doubling bifurca- 
tions around 0J2 and the resonance peaks around oj\. 
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Fig. 8. Left: largest nonzero Lyapunov exponents Ai^ as a function of v in the 
quasi-periodic doubling cascade. The first doubling points have been marked Dj, 
% = 1, 2, 3. Right: enlargement of the regions around Di and D2. Note that after the 
second doubling the Lyapunov exponents do not meet. 
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Fig. 9. Left: largest nonzero Lyapunov exponents Ai t 2 as a function of s for the 
model ODE (|15p with e = 0.001. The doubling points have been marked as in 
figure Right: enlargement of the region around T>3^- The bold lines denote the 
exponents for the perturbed model at e = 0.001. For comparison the exponents for 
the unperturbed case (e = 0) are shown with thin lines. Note that the exponents of 
the perturbed system do no meet between the doubling points. 
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